% ML TAU 2013 final project script
% Alt 6:
% claculate the hellinger norm for each pair of pictures in the data, 
% and learn the threshold seperating the results using SVM 

function gogo_alt6()


    %
    % Initialization
    %
    close all; clear; clc; tic;
 

	libsvm_path = './libsvm-3.17/matlab';

	% Add path to the libsvm
	addpath(libsvm_path);

   ker_degrees = 1:4;%1:4;

   ker_types = 1:2;%0:3;

   c_degrees = -3:8;%-3:7;
    
    if ~exist('dataforproject.mat','file')
        error('Data file not found in current directory')
    end;
    
    
    
    % Saving memory. Loading vars on a need to load basis only throughout
    load('dataforproject.mat','X1train','X2train','gidtrain','ytrain');
    fprintf('Training Data successfully loaded\n');
       
    
    X1train=X1train'; X2train=X2train';
    
    [m,n] = size(X1train);
    
    Data = zeros(m,1);

    for i=1:m
        
        P = X1train(i,:);
        Q = X2train(i,:);
        
        Data(i) = hellinger_norm(P,Q);
        
    end
    
    %[C D T Err] = svm_grid_search(Data,ker_types,ker_degrees,c_degrees);
    
    C=1.0000; D=1.0000; T=1.0000;
    Err = 1-do_svm(Data,T,C,D);
        
    fprintf('===\tResults\t===\n(C*=%f ; t*=%d ; deg*=%d)\n#Avg Acc: %f\n',C,T,D,1-Err);
    
    function [dist] = hellinger_norm(P,Q)
        
        r = 1/sqrt(2);
        
        dist = r*norm(sqrt(P/sum(P))-sqrt(Q/sum(Q)),2);
        
    end
    

function [accuracy] = do_svm(TX,t,c,d)
    
    accuracy = 0.0;

    for k=1:3

        % Training set
        X = TX((gidtrain~=k),:);

        y = ytrain(gidtrain~=k);

        % Create the model
        model = svmtrain(y,X,sprintf('-t %d -g  1.0 -c %f -d %d -h 0',t,c,d));

        % Testing set
        X = TX((gidtrain==k),:);
        y = ytrain(gidtrain==k);

        % Predict
        predictedY = svmpredict(y,X,model);

        results = sum(y.*predictedY > 0) / size(y,1); 
        
        accuracy = accuracy+results;
    end
    
    accuracy = accuracy/3;
end


function [C D T Err] = svm_grid_search(Data,ker_types,ker_degrees,c_degrees)
    
    s = size (c_degrees);
    
    c_vec = power(repmat(2,s),c_degrees);   
 
    Err = 1.0;
    
for  t=ker_types
   
    for d=ker_degrees
        
        for c=c_vec
            
            e = 1-do_svm(Data,t,c,d); 
            if e < Err
                C=c; D=d; T=t; Err=e;
            end
        end
    end
end
end

end


